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Abstract 



We present a novel density functional for liquid ^He, properly account- 
ing for the static response function and the phonon-roton dispersion in the 
uniform liquid. The functional is used to study both structural and dynam- 
ical properties of superfluid helium in various geometries. The equilibrium 
properties of the free surface, droplets and films at zero temperature are cal- 
culated. Our predictions agree closely to the results of ah initio Monte Carlo 
calculations, when available. The introduction of a phenomenological velocity 
dependent interaction, which accounts for backflow effects, is discussed. The 

spectrum of the elementary excitations of the free surface and films is studied. 
PACS number:67.40 
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I. INTRODUCTION 



The aim of the present work is to present a density functional theory, which treats 
static and dynamic properties of hquid hehum on the same ground, and is accurate at the 
microscopic level. 

The investigation of the properties of liquid helium in different geometries has a long 
story. Extensive work has been devoted to the search for and the understanding of su- 
perfluid effects in finite systems like helium droplets, or in quasi two dimensional systems 
like helium films on solid substrates, or within porous materials. Liquid helium exhibits 
very peculiar properties, such as the propagation of surface modes in the short wavelength 
regime, quantum evaporation of atoms produced by rotons impinging on the surface, wet- 
ting and prewetting transitions on solid substrates. The nucleation of bubbles at negative 
pressure and of quantized vortices near walls are further examples of interesting phenomena 
where both the inhomogeneity of the liquid and quantum correlations play an important 
role. These phenomena, among others, make liquid helium particularly appealing from the 
viewpoint of quantum many body theories. Theory and experiments, however, have not 
yet a satisfactory overlap. One difficulty comes from the fact that ab initio calculations are 
still hard to implement for inhomogeneous systems. On the other hand, phenomenological 
theories, which are quite successful in describing macroscopic properties, are not always 
adequate enough to investigate the behaviour of the system on the scale of interatomic 
distances. A major progress in this direction has been recently made in the framework of 
density functional theory (DFT). 

Several density functionals have been developed in the last years for applications to 
quantum fluids. The method consists of writing the energy of the many body system as 
a functional of the one-body density; from a given functional one extracts the equilibrium 
state, by minimizing the energy with respect to the density, as well as the excited states, 
by solving time-dependent equations of motion. An accurate phenomenological density 
functional for liquid ^He (Orsay-Paris functional) has been recently introduced it has 
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proven to be quite reliable in different contexts, such as helium films and wetting phenomena 
0-0, vortices in bulk hquid ||^ or droplets For dynamical properties it gives predictions 
which are close to the Feynman approximation for elementary excitations. In particular it 
does not account for backflow effects |0. 

In the present work we present a new functional. We follow the same ideas which lead to 
the Orsay-Paris functional, taking a similar two body interaction of Lennard- Jones type and 
including phenomenologically short range correlations. The main difference with respect 
to the Orsay-Paris functional is the addition of two new terms which allow for a better 
description of both static and dynamic properties on the scale of the interatomic spacing. 
These terms are: 

• a non local term depending on gradients of the density, which allows to reproduce the 
experimental static response function in a wide range of wave vectors, as well as its 
pressure dependence; 

• a term depending on local variations of the velocity field (backflow effects), which 
allows to reproduce the experimental phonon-roton dispersion in bulk. 

The velocity dependent term has a form similar to the one already introduced by Thouless 
T0| , who studied the flow of a dense superfiuid. The idea is to model backflow effects, which 



are important at small wavelengths, by a non-local kinetic energy term. 

With these ingredients the predictions for several properties of non uniform systems are 
significantly improved with respect to the ones given by previous functionals. The accuracy 
of the density functional theory is now comparable to the one of ab initio Monte Carlo 
calculations. The resulting approach, although phenomenological, represents a powerful 
and accurate tool whenever ab initio calculations become hard to implement. It allows one 
to investigate a wide variety of systems of different sizes and in different geometries, from 
few atoms to the bulk liquid, with limited numerical efforts. 

The paper is organized as follows: in Sect. |I| we introduce the density functional for 
static calculations. We emphasize and motivate the differences with respect to the Orsay- 



Paris functional. The results for the equilibrium properties of the free surface, droplets and 
films are given in Sect. In Sect. |IV| we discuss the application of the density functional 
theory to dynamics, showing the connection with the Feynman approximation and with the 
formalism of the Random Phase Approximation. We introduce a phenomenological current- 
current interaction which makes the density functional quantitative in the description of the 
spectrum of excited states. We discuss the role of this new term using sum rule arguments. 
Finally, in Sect. we present the results for the excited states of the free surface and of 
films. 



II. GROUND STATE CALCULATIONS 

In the density functional formalism at zero temperature, the energy of a Bose system is 
assumed to be a functional of a complex macroscopic wave function 

E = Jdrn[^,^*] . (1) 

The function \1/ can be written in the form 

vl/(r,t)=$(r,t)expQ^(r,t)) . (2) 

The real function $ is related to the diagonal one-body density by p = while the phase 
S fixes the velocity of the fluid through the relation v = {l/m)VS. In the calculation of 
the ground state, only states with zero velocity are considered, so that the energy is simply 
a functional of the one-body density p(r). A natural representation is given by 

E = JdvHolp] = Ec[p] + Jd^^i^Vpf , (3) 

where the second term on the r.h.s. is a quantum pressure, corresponding to the kinetic 
energy of a Bose gas of non uniform density. The quantity Ec[p] is a "correlation energy"; 
it incorporates the effects of dynamic correlations induced by the interaction. 

Ground state configurations are obtained by minimizing the energy of the system with 
respect to the density. This leads to the Hartree-type equation 
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{-|^V^ + f^[Ar]|vpW = /^VpW , (4) 
where U[p, r] = SEc/Sp{r) acts as a mean field, while the chemical potential p is introduced 
in order to ensure the proper normalization of the density to a fixed number of particles. 

For a weakly interacting Bose system the expression of Tio can be derived on a rigorous 
basis, yielding the well known Gross-Pitaevskii theory [|ri[] . Since liquid helium is a strongly 
correlated system, such a derivation, starting from first principles, is not available. One 
then resorts to approximate schemes for the correlation energy. Krotscheck and cowork- 
ers |T2|-[T^, for instance, have developed a variational hypernetted chain/Euler-Lagrange 
(HNC/EL) theory in which the one-body effective potential U is evaluated using the micro- 
scopic Hamiltonian. An alternative approach consists of writing a phenomenological expres- 
sion for the correlation energy, whose parameters are fixed to reproduce known properties 
of the bulk liquid. 

A simple functional was introduced in Ref. |^5[|l6l to investigate properties of the free 



surface and droplets of both ^He and ^He . The correlation energy, in analogy with the 



formalism of zero-range Skyrme interactions in nuclei [T^ , was written as 

E,[p] = Jdr 



-/ + Ip'"-' + d{Vpf 



(5) 

where 6, c and 7 are phenomenological parameters fixed to reproduce the ground state energy, 
density and compressibility of the homogeneous liquid at zero pressure, and d is adjusted 
to the surface tension of the liquid. The first two terms correspond to a local density 
approximation for the correlation energy, and non local effects are included through the 
gradient correction. The local density approximation is currently used to describe exchange- 
correlation energy in electron systems; its use for liquid helium does not provide quantitative 
results |18|. Non-local terms, like (Vp)^ in (|^), turn out to be crucial for the description 
of inhomogeneous liquid helium. Non-locality effects have been included in DFT in a more 
realistic way by Dupont-Roc et al. who generalized Eq. accounting for the finite 
range of the atom-atom interaction. The resulting functional has proven to be reliable in 
describing several inhomogeneous systems 



In the present work we follow the spirit of Ref. We use a similar two-body finite range 
interaction screened at short distances, and a weighted density (or "coarse grained" density) 
to account for short range correlations. The most important feature of this approach is that 
the static response function of the liquid can be reproduced even at finite wave vectors g, 
up to the roton region. The static response function fixes the linear response of the 
system to static density perturbations and is a key quantity in density functional theories. 
It is easily calculated from functional @ by taking the second functional derivative of the 
energy in g-space |]19| : 

A major advantage of liquid helium is that x(q') is known experimentally. It is related to 
the inverse energy weighted moment of the dynamic structure function S{q,uj), measured in 
neutron scattering, through the relation 

X{q) = -2m_i(g) , (7) 

the n-moment of S{q, u) being defined as 

mn{q) = / du S{q,uj){huoY . (8) 
Jo 

The experimental data for in the uniform liquid at zero pressure pO| are shown in Fig. |l| 
(circles). The prediction of the zero-range functional of Ref. [|l5l (i.e., with the correlation 
energy given in Eq. (^) and of the finite range Orsay-Paris functional [|l| are also shown as 
dotted and dashed lines respectively. The g = limit is fixed by the compressibility of the 
system, which in both cases is an input of the theory. It ensures the correct behaviour in the 
long wavelength limit and, consequently, the correct description of systems characterized by 
smooth density variations as happens, for example, in the free surface or in helium droplets. 
The height of the peak of the static response function in the roton region, g ^ 2 A~^, is 
important in characterizing structural properties on the interatomic length scale (for instance 
the layered structure of helium films). Therefore, the first important idea is to improve on 
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the Orsay-Paris functional in order to better reproduce the experimental peak of in the 
roton region. 

In the present work, the function Tio, entering functional is taken in the form 



H„= |^(VV7;)' + i/<ir'p(r)V,(|r-r'|)p(r') + |p(r)fe)2 + |p(r)(ft)» 

The first term on the r.h.s. is the quantum pressure, as in Eq. (^. The second term 
contains a two-body interaction VJ, which is the Lennard- Jones interatomic potential, with 
the standard value of the hard core radius of 2.556 A and of the well-depth 10.22 K, screened 
at short distance (V/ = for r < /i, with h = 2.1903A). The weighted density p is the average 
of p(r) over a sphere of radius h: 

p(r) = |rfrn,(|r-r'|)p(r') , (10) 

where Ilh{r) = 3(47r/i'^)^^ for r < h and zero elsewhere. The two terms containing p, 
with the parameters C2 = -2.411857 x 10^ K and Cg = 1.858496 x 10^ K A^ account 
phenomenologically for short range correlations. All these terms have a form similar to 
the Orsay-Paris functional ||l|. The Lennard- Jones potential is here screened in a simpler 
way, avoiding the fourth power for r < h, and the dependence on p is slightly different. 
The effects of these minor changes will be discussed later. The last term in Eq. (|^) is 
a completely new term; it depends on the gradient of the density at different points and 
corresponds to a non-local correction to the kinetic energy. The function F is a simple 
gaussian F(r) = 7r~^/^£~'^ exp(— r^/£^) with £ = 1 A, while = 54.31 A^. The parameters 
are fixed to reproduce the peak of the static response function in the bulk liquid. The latter 
can be easily calculated by using Eq. (^. One finds 

- x-\q) = ^ + pViiq) + C2(2n,,(g) + n^(g))/ 

+ 2c3(fl.(g) + Il^(g))p^ -^a.p(l-^)Vexpf-^) , (11) 

2m pos V 4 / 
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where Vi{q) and tlh{(l) are the Fourier transforms of the screened Lennard- Jones potential 
and the weighting function U^, respectively, while p is the bulk density. The resulting curve 
for the liquid at zero pressure [p = po = 0.021836 A~^) is shown in Fig. (solid line). The 
factor (1— p/ Pos), with pos = 0.04 A~^, is included in order to obtain a pressure dependence of 
the static response function close to the one predicted by Diffusion Monte Carlo simulations 
21| . For instance in the liquid near solidification (p = 0.02622 A^'^) one finds a peak of 



about 10 % higher than at zero pressure, and displaced by 0.1 to larger wavelengths. 
Finally, the quantity p(r) is again a weighted density, calculated using the gradient-gradient 
interaction function F as a weighting function: 

p(r) = Jdr'F{\r-r'\)p{r') . (12) 

Actually the density p(r) is very close to the particle density p(r) except in strongly inho- 
mogeneous situations (like helium adsorbed on a graphite substrate). For this reason one 
can safely replace p(r) with p(r) for the investigation of the free surface, helium droplets and 
films on weak binding substrates. Stronger constraints on the form of F should be provided 
by the study of the liquid-solid phase transition. 

In a uniform liquid of constant density the energy per particle, from functional (|,^, 
reduces to the power law 

where b = —718.99 KA'^ is the integral of the screened Lennard- Jones potential VJ. Pressure 
and compressibility can be derived directly by taking the first and second derivative of the 
energy: 

9 d E , dP 

The experimental values of the density, energy per particle and compressibility of the uniform 
system at zero pressure are used as input to fix the parameters h, C2 and C3. The resulting 
equation of state and the sound velocity, = — [mx(0)]~^, are shown in Figs. ^ and |^ 
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respectively. The comparison with the results of Monte Carlo simulations and with 



experimental data p3| , |2^ shows that the present density functional theory describes correctly 
the ground state of the bulk liquid ^He at all pressures. Comments about the equation 
of state in the limit of a quasi two-dimensional liquid will be given when discussing the 
structural properties of films. 

To conclude this section we emphasize again the main idea. The density functional (||,^ 
incorporates the correct long range behaviour of the interatomic potential, and accounts 
for short range correlations in a phenomenological way. The ground state of the uniform 
liquid is well reproduced at all pressures and, furthermore, the response to small static 
density perturbations, up to the roton wave length, is also correctly reproduced. Significant 
differences with respect to previous functionals are expected in the predictions of properties 
which depend on the behaviour of the fluid on the interatomic length scale. Some interesting 
examples are discussed in the following section. 



III. RESULTS FOR THE EQUILIBRIUM CONFIGURATION OF FREE 
SURFACE, DROPLETS AND FILMS 

A. Free surface 

To compute the density profile and the energy of a planar free surface of "^He at zero 
temperature one has to solve the integro-differential Hartree equation with the mean 
field U extracted from (^). Both p and U depend only on the coordinate orthogonal to the 
surface, so that the equation is one- dimensional. The integrals on parallel coordinates in 
the non-local terms of the functional can be written analytically. The numerical solution of 
non-linear Eq. (^) is obtained by means of an iterative procedure and provides the density 
p{z), from which the surface tension a can be calculated through 

a = j dz {Hq[p] - pp} . (15) 
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The bulk density is kept fixed to the experimental value at zero pressure, po = 0.021836 
A~'^. The density profile is shown in Fig. ^, together with the results given by the Skyrme 
functional of Eq. (|^) and the Orsay-Paris one. The surface tension is practically the same 
in the three cases. The present functional gives a = 0.272 K A~^, to be compared with the 
experimental values quoted in the literature: 0.275 K and 0.257 K A^^ [^]. The 

10%-90% surface thickness is approximately 6 A. The value 7.6 A, extracted from X-ray 



scattering data is slightly larger, being closer to the result of the zero-range functional 
of Eq. However, the experimental value is expected to depend on the form of the density 
profile used to fit the measured refiectivity. 

Notice that the density p{z) resulting from expression (|]) is not as smooth as the ones 
given by previous calculations. It exhibits small oscillations which appear as shoulders on 
the surface profile and asymptotically die in the bulk region. Oscillations of this type were 
predicted long time ago by Regge ^8[. In his theory the surface was treated as a source 
of elementary excitations producing static ripples on the density profile. The form of the 
ripples is connected to the behaviour of the static response function. The pronounced peak 
at the roton wavelength is associated with the tendency of atoms to localize in "soft sphere 
close packing". This tendency is opposed by the zero point motion of the surface, whose 
thickness is larger than the interparticle distance, and by compressibility effects, so that the 
density oscillations turn out to be quite small. The role of the static response function in 
characterizing the form of the surface profile is clearly seen on Figs. |I] and ^. Of course one 
expects that functionals describing the short wavelength behaviour of xW) in a correct way 
will give rise to more reliable predictions for the properties of the fiuid in the microscopic 
scale. This will be explicitly confirmed in the case of helium droplets, where our predictions 
can be compared with Monte Carlo simulations. 
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B. Droplets 



The solution of the Hartree equation in spherical symmetry with a fixed number of 
particles N provides the ground state of ^He droplets. Since the density functional approach 
is not time-consuming, it allows to compute the energy and the density profile of droplets 
in a wide range of sizes. The density p(r) for droplets with 8 < < 60 is shown in Fig. |^. 
Again one finds ripples on the surface profile; they are more pronounced than in the case of 
the free surface, since the droplets have a size of the order of a few interatomic distances and 
the "soft sphere close packing" tends to produce shell structures. This effect was already 
suggested by Rasetti and Regge ||29|, but subsequent theoretical calculations [|l^ ,^0H32| did 
not predict any clear and systematic oscillation in p(r). Only recently, sizable oscillations 
have been found in the Diffusion Monte Carlo (DMC) calculations by Chin and Krotscheck 
||33| . An example is given in Fig. |^, where we show the density profile of a droplet with 
70 particles. The solid line is the result of the present density functional calculation, while 



the DMC results of Ref. are represented by triangles. The DMC data exhibit more 
pronounced oscillations, but the presence of long-lived metastable states, slowing down the 
convergence in the Monte Carlo algorithm, cannot be completely ruled out [^]. Recently the 
same authors have found oscillations in p(r) even with a variational approach based on the 
HNC approximation (dashed line). Even though the HNC method underestimates the 
central density, it predicts oscillations with amplitude and phase in remarkable agreement 
with the ones of density functional theory. An even better agreement is found in the most 
recent DMC calculations by Barnett and Whaley ||3^ (circles), where the statistical error is 
significantly reduced with respect to Ref. fS^ . 

The detailed structure of p(r), though interesting from a theoretical viewpoint, is not yet 
observable experimentally with enough accuracy to distinguish between a smooth profile and 
a profile with small oscillations. It is thus important to explore the effects of the "soft sphere 
close packing" on the energy systematics, since the latter is related to the mass distribution 
of droplets in the experimental beams In Fig. ^ the energy per particle is given as 
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a function of A^. The results of the present work (full line) are compared with the ones 
of previous functionals (dashed and dotted lines), as well as with Monte Carlo simulations 
(dots from Ref. I^S] and crosses from Ref. |0). First we note that the accuracy of the 



density functional theory, compared with "ab initio" simulations, increases progressively, 
following the improvement in the prediction of the static response function x(g) in the 
microscopic region. The agreement between the results of our functional and Monte Carlo 
data is excellent. Second, as the energy is a smooth function of A^, helium clusters behave 
essentially as liquid droplets. Indeed the energy can be easily fitted with a liquid drop 
formula: 

^ = a^ + a,A^-i/3 ^ a^jV"2/3 + aoN~^ , (16) 

where the volume coefficient is the chemical potential in bulk, the surface energy is 
fixed by the surface tension, while Oc and oq can be taken as fitting parameters. The energy 
calculated with the density functional differs from the liquid drop fit by less than 0.02 K for 
all droplets with > 30. This seems to rule out apparently any shell effect in the energy 
systematics. However, the relevant quantity to investigate in this context is the evaporation 
energy [E{N — 1) — E{N)]. The latter is not as smooth as the energy per particle. Figure 
H shows the evaporation energy predicted by the density functional (solid line) and the one 
obtained with the liquid drop formula ([TgD (dashed line). The difference between the two 
curves is also shown (circles). When the difference is positive the droplets are more stable 
than it is predicted by the liquid drop formula. We note clear oscillations, having decreasing 
amplitude and increasing periodicity as a function of A^. The same kind of oscillations appear 
in the central density of the droplets, as seen in Fig. |[ Since the distance between two crests 
of the surface ripples is practically constant and the droplet radius goes approximately like 
N^^^, the period of oscillations of the central density as a function of A^, as well as the one 
of the evaporation energy, increases as N^^^. The predicted deviations of the evaporation 
energy from the liquid drop behaviour are rather small (less than 0.1 K). Unfortunately this 
value is smaller than the temperature of droplets in available experimental beams (about 
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0.4 K). 



C. Layering and prewetting transitions in films. 

In order to find the equilibrium state of liquid helium on a solid substrate, we add the 
external helium-substrate potential KMfe(r) to the mean field U in the Hartree equation (^. 
We assume the substrate to be flat, avoiding the problem of possible corrugations. This 
approximation is certainly valid for weak-binding substrates, such as the alkali metals, for 
which the cloud of delocalized electrons is expected to smooth out the potential along the 
substrate plane. In this case the Hartree equation is again one dimensional. The substrate- 
helium potential is taken here as a 9-3 potential of the form 

where C3 is the Hamaker constant and D the well depth of the potential. The values taken 
from for various substrates are gathered in Table |. 



The energetics of helium films on various surfaces indicate that two types of phase transi- 
tion in film growth can take place, depending on the strength of the substrate potential. On 
strong and medium binding surfaces, which create large local pressures (one or two layers 
may become solid) the growth of the first liquid layers proceeds via layering transitions, first 
described for helium on graphite in Ref. |3^. When considering weak binding substrates, 
however, these layering transitions are no longer present. One enters a new regime of film 
adsorption where prewetting transitions take place at T = 0, as analyzed in details in Ref. 0. 
By convention, the term layering transition is reserved usually to the case of a first order 
transition in film thickness involving one layer only. The occurence of such transition is 
not linked to the question of wetting; it is essentially related to the nature of the quasi-2D 
system. To the contrary, a prewetting transition involves a jump in film thickness which can 
take any value. The notion is intimately connected to that of wetting, since a prewetting 
transition is the continuation of the wetting transition off coexistence. 
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The various cases are best illustrated by considering, for a given substrate, the de- 
pendence of the chemical potential fi on coverage, as shown in Fig. |^. A negative slope 
{dfi/dN < 0) indicates an unstable range of film thickness, and the transitions are deter- 
mined by Maxwell constructions. The figure shows the results for a H2 substrate and for the 
alkali metals. For H2, we have used a parametrization of the helium-substrate interaction 
proposed in Ref. ||40[| : 



, , 900000 15000 435 , , 

VH,iz) = —, z^--r ' (18) 



Z5 



where z is in A and in K. This form is fitted to the results of Pierre et al. |^T| and gives 



a well depth of 33 K. The binding energy of one ^He atom is 16.4 K, in good agreement 
with the experimental determination of Paine and Siedel [0, which is 16±2 K. For this 
potential, two layering transitions occur. The mechanism by which these transitions are 
produced can be summarized as follows: on a strong substrate, liquid helium forms well 
defined layers which are approximately independent quasi-2D systems. As 2D helium is a 
liquid (Monte-Carlo simulations |^ give a binding energy of 0.8 K at an equilibrium density 
of 0.043 A~^) the formation of each of the first layers exhibits a quasi-2D condensation. It is 
important to realize that the two aspects are necessary : on one hand the layering transitions 
would not occur if quasi-2D ''He were a gas; on the other hand we shall see below that they 
also disappear if the substrate is not strong enough to produce a sufficient layering of the 
fluid, as is the case with the alkalis. Fig. ^ represents the growth of a helium film on H2 
characterized by two regions of instability or metastability. 

On sees from Fig. ^ that the Orsay-Paris functional produces rather smooth curves and 
in particular misses the layering transitions besides the first one. This fact is due to two 
deficiencies of the model, related to the mechanism described above, namely i) the peak of 
the density-density response function is underestimated, as we have seen, by a factor of 
almost 2 and ii) the binding energy of the quasi-2D system is too small. Both defficiencies 
are corrected with the new functional. 

Let us now turn to the results for the alkalis. Although fi{N) has still some structure, 
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Mg appears as a limiting case where the layering transitions tend to disappear, except for 
the first one. Interestingly, one sees that with decreasing strength of substrate potential, 
it is this first transition which becomes larger and larger in amplitude and thus turns into 
the prewetting transition. The physics here is no longer that of the quasi-2D system, but 
that of wetting and prewetting : for a given substrate, the thinnest stable film is such 
that the energy cost of forming two interfaces -one with the substrate, one free surface- is 
compensated by the energy gain of placing the fluid in the attractive potential of the wall. 
The limit between wetted and nonwetted substrates is obtained when stability is obtained 
only for an infinitely thick film. Notice that the predictions of the Orsay-Paris functional 
and of the present one become similar for these weak binding surfaces. Using the original 
values of the substrate potential parameters, one still finds that the three alkalis Cs, Rb and 
K are not wetted. The new functional slightly favors wetting with respect to the Orsay-Paris 
one. For example, the contact angle calculated for Cs is reduced by 3 degrees; also, on a Na 
substrate, the prewetting jump is reduced to 3.3 layers, compared to 5.2, and metastable 
films are found to exist down to 1.2 layers. 

Layering growth can be also seen in the adsorption isotherms. In Fig. ^ is plotted the 
isotherm T = 0.639 K for helium on graphite. For comparison with experimental data |^ 
we use the ideal gas formula 



P 

— = exp 



(19) 



where /iq = —7.15 K is the bulk chemical potential, Pq the saturating vapor pressure at T. 
The graphite substrate is a test for the model in a highly inhomogeneous situation. The two 
first layers are known to be solid. To first approximation, the effect of localization of the 
atoms in the plane parallel to the substrate can be ignored. Hence, we have treated the two 
solid layers as the liquid. The experimental data show a clear staircase structure, associated 
with the progressive filling of layers. The results of the present functional (solid line) exhibit 
a similar pattern, with steps of amplitude and phase close to the experimental ones. 
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IV. DYNAMICS 



In Sect. we wrote the Hartree equation (Q) for the ground state of the fluid, which 
corresponds to the minimization of the energy E in Eq. (|I]) with respect to the density. 
This formahsm is generahzed to the study of dynamical properties using the least action 
principle: 



rt2 

5 J dt I dr 



. (20) 



The equations of motion for the excited states of the fluid can be derived by making variations 
with respect to \1' or \E'*. One finds a Schrodinger-like equation of the form 

(/f-/i)vi/ = ^^^xi/ , (21) 

where H = 5E/5'^* is an effective Hamiltonian. We linearize the equation by writing 

^(r,t) = ^o(r) + 5^(r,t) (22) 

where \I/o(r) refers to the ground state. The Hamiltonian H then takes the form 

H = Ho + 6H . (23) 



The static Hamiltonian 



Ho = -^V' + U[p,r] (24) 



which appeared already in Eq. (H) determines the equilibrium state \E'o(r) = \J p{'^)- The 
term 5H is linear in 5'^ and accounts for changes in the Hamiltonian induced by the collective 
motion of the system. Since H depends explicitly on the wave function the Schrodinger 



equation (^) has to be solved using a self-consistent procedure, even in the linear limit 
considered in the present work. 

The formalism here described corresponds to a time dependent density functional 
(TDDF) theory which, in the linear limit (P^D, coincides with the Random Phase Approxi- 
mation (RPA) for a Bose system. This theory, which is basically a mean field theory coupling 
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one-particle/one- hole configurations, is suitable for describing collective (one-phonon) states, 
but not multiphonons excitations. 

A completely equivalent formulation of the equations of motion can be obtained by using 
the canonically conjugate variables p and 5, defined in Eq. (|^). The least action principle 
takes the form 



rt2 

6 dt I dr 



dS 

T-Co[p] + K[p, v] - /xp + p— 



, (25) 



where we have separated the velocity dependent part of the functional from the velocity 
independent one. After variations with respect to p and S, one finds two coupled equations 
of the form 

% + -^r^fdv{n.[p,^r]} = (26) 
ot m ov J 

^ + j-p Jdr {Holp] + n,[p, v] - /ip} = , (27) 

which can be viewed as a generalization of the equation of continuity and the Euler equation. 
In the hydrodynamic limit, Ti^, is given by the usual kinetic energy term 

K[p,v] = ^|vp , (28) 

so that ( P^ ) leads to the usual conserved current 

Jo(r) = p(r)v(r) = ^VS(r) . (29) 

The present formalism, with Ti^ given by (^) and with different choices for TYq, has 
been already applied to helium droplets and films However, there are two major 

shortcomings which make the results of those calculations not satisfactory from a quanti- 
tative point of view. First, the static response function is not well reproduced in the roton 
region. Second, the theory does not account for backfiow effects. 

To understand this point better, let us discuss the results of the TDDF in bulk. In 
this case the density is constant and (5\E' can be expanded in plane waves, corresponding to 
the propagation of phonon-roton excitations. In the absence of backfiow effects, the only 
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dependence of the energy on the velocity field comes from the hydrodynamic limit of 7i„ 
(see Eq. (pSf ) and the resulting dispersion relation for the phonon-roton mode takes the form 

Result (|30|) can be also written in terms of the moments (H) of dynamic structure function 



S{q,uj). Indeed, one has 



M,) = f^V" . (31) 



where the energy-weighted moment 

|"oo h^q'^ 

^lio) = I duj S(q,uj)huj = (32) 

JO 2m 

coincides with the well known f-sum rule, while the inverse energy- weighted moment m_i is 
related to the static response function by the compressibility sum rule (^. 

Density functionals having the form ( [28| ) for Ti.^ exactly fullfil the f-sum rule. Thus, 
differences in the predictions for the dispersion uj{q) come only from the quantity x(Q'); 
which is fixed by the static part 7io[p] of the functional. In Fig. ^ we show the results 
for the phonon-roton dispersion obtained with the Orsay-Paris functional (dashed line) and 
functional (§) (dot-dashed line). Both curves overestimate significantly the experimental 
phonon-roton energy (points). The difference between the predictions of the two functionals 
is clearly understood by looking at Fig. |l| and Eq. (0): starting from the hamiltonian 
density of Eq. (P), one obtains the full peak of x(q') at the roton wavelength, thus predicting 
a lower roton energy. The remaining gap between theory and experiment is mainly due to 
the role of multiphonon excitations in the f-sum rule. In fact, the analysis of the spectra 
of neutron scattering experiments |2y] shows that the collective mode gives only a fraction 
(^ 1/3) to the f-sum rule in the roton region, the remaining part being exhausted by high 
energy multiphonon excitations. On the contrary, due to the uj~^ factor in the integrand, the 
collective mode almost exhausts the compressibility sum rule for wavelengths up to about 2.2 
A~^. This means that in order to have an accurate prediction for uj{q), only the single mode 



contribution to the f-sum rule should be included in the moment mi{q) entering Eq. (^TJ). 
Therefore, one concludes that, once the static response function is properly accounted 



for, the dispersion law (^) given by TDDF theory provides an upper bound to the exact 
energy of the phonon-roton mode. 

It is instructive to compare the above predictions with the results of the so called Feyn- 



man approximation for collective excitations , often used in dynamic calculations for non 



uniform ^He states |ll3| , |32| , |33| , ^ . In the bulk the Feynman dispersion law takes the form 



where the static structure factor is related to the non energy-weighted moment of S{q, u) 
through the equation 

S{q) = mo{q) = JdujS{q,uj) . (34) 
General properties of the moments m„ permit to prove the following inequality: 



< ^ (35) 

holding at zero temperature. This implies that the dispersion law given by Eqs. ( PDIPH ) 
provides an upper bound closer to the exact dispersion law with respect to the Feynman 
approximation. This is explicitly shown in Fig. |12|. 

Concerning the static structure factor S{q), one should note that it can not be properly 
accounted for by the TDDF theory developed above. The reason is that, in the bulk liquid, 
this theory represents a single mode approximation with the dynamic structure function in 
the form 

5^^(g, u) = -^^6{uj - udf) (36) 
ZmujDp 

and the dispersion law udf given by Eq. (pT]). While Eq. ( p6D reproduces exactly the sum 
rules mi and m_i, it yields the approximate expression 



5^^(g) = ^m,{q)m^,{q) = sj 



Am 
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^ \x{q)\ (37) 



for the static structure factor. Actually result (|37|) provides an upper bound to the exact 
value of S{q). The difference is again due to the role of multiphonon excitations. 

One of the purposes of this work is to make the TDDF theory more quantitative. The 
idea is to realize that the expression of Ti^ given by Eq. (^8|) comes from a many-body wave 
function of the form 



^Ar(ri, ...tn) = exp 



|^^(ri,...r,v)| , (38) 



where s is real. Clearly, this wave function does not take into account short range phase 
correlation, and describes correctly the superfluid motion in the hydrodynamic limit only. 
The exact wave function should be expressed with the more general phase: 

s(ri,r2,...rN) = + Zl'^2(ri,rj) + ... + ^ s^v (r^, r^, r^, ...) . (39) 

i i<j i<j<k... 

The average kinetic energy obtained from this wave function has non-local contributions 
coming from S2, ■■■Sn. Now, if we make the assumption that those terms can be expressed 
only with the two canonically conjugate variables p and S, then we are led to add a non-local 
velocity dependent term to 7^t,[p, v]. Indeed, this procedure has been proposed a long time 
ago by Thouless [|TU| in the study of the flow of a dense superfluid, but in his article, helium 



was treated as an incompressible liquid. In the present work, we incorpore the suggestion 
of Thouless by introducing the most general quadratic form 

= -p(r)|v(r)|2 + Jdr.dr^dr, [v(r) - v(ri)]G(p; r, n, r2, r3)[v(r2) - v(r3)] , (40) 

where we have explicitly exhibited Galilean invariance. Generally speaking, G is a tensor 
which takes non zero values on scales |rj — rj \ of the order of the inter-particle distance. So 
far no microscopic derivation has been found for this functional form. Our phenomenological 
approach consists in keeping only diagonal terms as follows: 

n. = '^p{r)\v{r)\'-jJdr'Vj{\r-r'\)p{r)p{r') [v(r) - v(r')f , (41) 

and fixing the effective current-current interaction Vj to reproduce known properties in bulk. 
The new term plays the role of a non local kinetic energy. Now, the conserved current is no 
longer Jq; rather, Eq. (p6D leads to the current 
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J(r)=Jo(r)+JB(r), (42) 
JB(r) = p(r) J p(r') [v(r) - v(r')] Vji\r - v'Ddv' (43) 

Physically, Js(r) acts as a backfiow which depends on the velocity and the density in the 
vicinity of point r. As expected, its contribution vanishes when many-body phenomena are 
not present (p ^ 0). 

This form of 7i„ lowers the value of the energy weighted moment rrii (g) predicted by the 
density functional theory. In bulk one now finds 

m,{q) = ^{l-p [VjiO) - Vjiq)] } , (44) 

where Vj(g) is the Fourier transform of the current- current interaction Vj(r). Notice that 
the expression for the static response function x(g) does not change, since it is entirely fixed 
by Tio- The phonon-roton dispersion in bulk is still given by the ratio (pi]), so that the 
dispersion law is given by 

This relation can be used to fix Vj(r) in order to reproduce phenomenologically the experi- 
mental phonon-roton dispersion. We have chosen the simple parametrization 

= (711 + 7i2?^^) exp(-air^) + (721 + 722?^^) exp{-a2r^) , (46) 

where the parameters are given in Table The corresponding dispersion relation is shown 
in Fig. |l^ (solid line). The roton minimum is at go = 1-92 and the roton energy is 
A = 8.7 K. The pressure dependence of the dispersion relation turns out to be also well 
reproduced. At P = 24 atm, for instance, the roton minimum is displaced at g = 2.01 
A~^ and the roton energy is A = 7.4 K, close to the experimental values g = 2.05 A~^ and 



A = 7.3 K pp 



The energy weighted moment mi{q) is shown in Fig. O, where the density functional 



result ( P^ ) is compared with the experimental data for the collective contribution to the f- 
sum rule . The fact that the new functional no longer satisfies the f-sum rule (|32D points 
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out in a clear way that the TDDF theory does not account for multiphonon excitations. 
The new current-current term in the functional ( |^ ) changes also the expression for the 



static structure factor, which is still given by S^^{q) = ymi(g)m_i(g) as in Eq. (p7D, but 
with the new mi{q) moment (H). The resulting S^^{q) is no more an upper bound to the 



exact S{q)] conversely it turns out to be close to the experimental one-phonon contribution 
to the total S{q). A similar separation between collective and multiphonon excitations, in 



the context of linear response theory, has been developed by Pines and collaborators ^2 
Finally we note that the dispersion law (|45|) is not affected by the new term in the g — 
hydrodynamic regime, where it gives the usual sound velocity = —[mx{0)]^^. This is an 
important feature ensured by Galilean invariance. 

In conclusion, the complete Orsay-Trento functional has the form 

E = Jdr {no[p]+n.[p,^r]} , (47) 

where TCq and Hy are given in Eqs. (|^) and (^) respectively. The term vanishes in 
the ground state calculations of the previous sections. On the contrary, it is crucial in the 
calculation of the dynamics. Both the experimental static response function and the phonon- 
roton dispersion in bulk are taken as input to parametrize the functional. The theory can 
then be applied to study the excited states of inhomogeneous systems. This is the purpose 
of the next section. 



V. EXCITED STATES OF THE FREE SURFACE AND OF FILMS 

A. Equations of motion. 

A detailed study of excitations using the Orsay-Paris functional can be found in Ref . . 
The same method of resolution is followed here. The systems under consideration have 
translational invariance parallel to the x — y plane, hence the ground state wave function 
depends only on one coordinate, which is taken as the z-coordinate orthogonal to the free 
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surface or to the plane of a film. The excited state wave functions can be expanded in plane 
waves as follows : 



where ry denotes from now on a two-dimensional vector parallel to the surface of the liquid. 
The quantities ak,fe are small amplitudes, the factor 1 / V AL is a normalization constant {A 
is the area of the sample, L is an arbitrary length) and b labels the various branches of 
excited states. The functions $i(-s) and $2 (-2), without any loss of generality, can be chosen 
real. This formalism has been already used in Refs. P,p!8|. 
It is convenient to introduce the quantities 



(49) 



Then the phase S can be written as 

S{r\\,z,t) = 



ih ^ 1 

2VAL^, ^0 



and the density takes the form 



(50) 



p{^,z,t) = po{z) +5p{rii,z,t) 



(51) 



where po = I^E^ol is the ground state density, while 

5pi^,^,t) = -l=^vl/o(z){a^,<|.+,(z)e-^(-'=>^*-''-r|^ + 



c.c 



The equations of motion (pB]) and (^) assume the form 



(52) 



(53) 



1^k<^k,b + ^{%b} = ^^k,b<^tb 
The eigenvalues huk^b satisfying Eq. ( p3D are the energies of the collective modes of the 
system. The symbol Vk denotes a differential operator of second order, while W and U are 
integral operators. They are defined by 
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Vk = -- — + — h (7[po, z\ + Vsub[z) - H 



2m dz"^ 2m 



6p{r)5p{r') 



JMr\\-r\\') 



dr 



.»k-(r||-r||') 



(54) 
(55) 

(56) 



_SSir)SS{v') 

where Ey is a functional containing only the non-local part of Ti^. A simple way of solving 
Eq. ( [53| ) is to expand it on the basis of eigenstates of the static one-body Hamiltonian. This 
leads to a matrix equation that is easily solved numerically by a direct diagonalization. 

A quantum hydrodynamical formalism is also obtained by making the substitution 
{— iak.b} Ok,fe, where ak,6 is the operator of creation of one phonon which statisfies the 
Bose commutation rules 



(57) 



From the completeness relation 



(58) 



one recovers the equalities 



(59) 



[p(r), Sir')] = zM(r - r') [p(r), p(r')] = [5(r), 5(r')] = , 

Then, after a quadratic expansion of the total energy and using the properties of Eq. (p3[), 
one finds the Hamiltonian 



(60) 



where Eq is the energy of the ground state. This diagonalization is obtained provided the 



$^^'s are normalized according to 



(61) 



where the orthogonality appears as a consequence of (|53|). Higher order expansions of H 
in terms of density and phase fluctuations would give rise to interactions between quasi- 
particles. 
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The transition density associated with the solution ([5^ ) is given by (the vector q denotes 
(k,g.)) 

K)m = "'/l / ^o(^)$m(^) exp{tq,z)dz (62) 
and the dynamic structure function can be evaluated by means of the definition 



b 



^(cj-cuk.b) • (63) 



In order to compare the theoretical strength with the experimental results, we have intro- 
duced a width of the order of the experimental resolution, by substituting the 5-function of 
Eq. ( |63D with a normalized gaussian, as done in Ref. ||53[ . 



B. Dynamics of the free surface. 

We have done calculations in slab geometry, i.e., liquid between two parallel surfaces 
at a distance L much larger than the surface thickness (typically 50 to 100 A). The slab 
geometry is a good approximation to a semi-infinite medium for wave vectors larger than 
1/L. The dispersion relation of the various modes, extrapolated to L ^ oo, are shown in 
Fig. [T^. The lowest branch corresponds to a wave function $+(2;) localized in the surface 
region. In the long wavelength limit, it coincides with the hydrodynamical surface wave 
called ripplon, whose dispersion is 

^^^k) = ^ , (64) 

where a and poo are the surface tension and the bulk density respectively. For wave vectors of 
the order of 0.5 A~^, the dispersion relation starts deviating from the hydrodynamical limit, 
until its curvature eventually changes sign, due to a coupling with bulk modes. Our curve 
reaches the value of the roton energy A = 8.7 K at about k = 1.15 A~^. A first experimental 
evidence of a similar deviation from the hydrodjTiamic law came from measurements of the 
surface entropy More recently the dispersion of surface modes has been measured in 



neutron scattering experiments on helium films [g^. The experimental data are shown in 
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Fig. [T^ as squares, the error bars accounting for the spreading of the data for different 
coverages (3.5 to 5 layers of hehum on graphite). The agreement with the calculated values 
is satisfactory. 

Below A the surface modes are undamped, while above A they couple with the continuum 
of bulk modes (rotons with negative and positive group velocity) propagating at different 
angles {qz 7^ 0). This results in a spreading of the strength associated with the surface 



modes. Actually the spreading predicted by our theory is small. In Fig. ^ we show the 
dynamic structure function for scattering at grazing angle (g^ = 0) on a slab of thickness 
L = 50 A. The strength of the lowest surface mode is well localized not only below but also 
above A, even though it is partially distributed among bulk modes coupled to ripplons. The 
position of the peak of the ripplon mode above A is shown also in Fig. |1^ as a dashed line. 

We obtain also a second branch of surface excitations, lying in between the ripplon and 
the phonon-roton modes. It is visible in Fig. |15| as a small bump, which is close to the 
phonon peak for k ~ 0.3 A^^ and stays almost parallel to the ripplon dispersion for larger 
k. The relative strength of this mode is larger in thin films; in that case, the experiments 
seem also to support the existence of such a surface mode. 

A more detailed analysis of the results of the present density functional theory for the 
excitations of the free surface is given in Ref. |5^; in that work, a general discussion of the 
mechanism of hybridization between ripplons and rotons is presented, and some properties 
connected with reflection and evaporation of bulk excitations impinging on the surface are 
also discussed. 



C. Dynamics of films. 

In order to illustrate the results obtained with the present model, we have chosen the 
case of a II2 substrate, since it has been the subject of both experimental and theoretical 
investigations PJ^. Besides finite size effects, the interesting features to be expected are 
linked to the layering of the liquid near the substrate. Submonolayer superfluidity has been 
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observed on H2 [0], indicating that helium remains hquid close to the substrate. However 



the well-depth of the helium-hydrogen potential is rather large (33 K), which produces a 
well defined layering of the liquid. 

In a thin film, the long wavelength limit of the lowest surface mode is no longer a ripplon, 
for which the restoring force originates from the surface tension, but rather a third sound 
excitation. The restoring force is here given by the substrate potential. The third sound 
speed is obtained as the hydrodynamic limit (k — ^ 0) of Eq. (^3|): 

-L-iV^, (65) 

where N is the coverage of the film. Indeed, one can verify that for small momentum, 
$^ ~ \l>o, so that the contribution due to U vanishes. With increasing film thickness, since 
the chemical potential varies as so does c|^. 

The results for c^s (Fig. |16D, show strong oscillations as a function of coverage. This 
is due to the layered structure of the film, which is refiected into the coverage dependence 
of the chemical potential (see Sect. [Ill (J| above). Films for which a uniform coverage is 
unstable have a constant chemical potential characteristic of a layering transition, the range 
of which is determined by a Maxwell construction as discussed above. This causes the third 
sound velocity to drop dramatically to zero. The structure is more marked than with the 
Orsay-Paris functional 0, which missed the layering transitions besides the first one (see 
Sect. [Ill Q) . Experimentally, there is up to now no indication of the layering transitions. The 



third sound velocity does show oscillations [57|, however less marked than calculated here. 



although the measurements were done at a temperature of 0.18 K, which is expected to be 
lower than the critical temperature of the layering transitions. It is also possible that surface 
inhomogeneities of the H2 surface are able to smooth out the dips seen in the calculations 
of for a perfectly fiat surface. 

Two typical spectra are shown in Fig. |1^ for two values of helium coverage (A^ = 0.24 
and N = 0.48 A~^) on II2 substrate. For wave vectors in the range 1.5 to 2 A~^, the low 
lying excited states are modes confined in the first layers of the fiuid close to the substrate. 
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Their nature has been discussed in details in Refs. PJ53[]. Their formation reflects the fact 
that the substrate-hquid interface is nearly solidified, so that low energy excitations appear, 
trapped in each of the first layers, with a momentum corresponding to the interparticle 
distance. The minimum 2D roton energy was found to be about 12 K with the Orsay-Paris 
functional. It reduces to 5.6 K with the present one. The combined effect of the evolution 
of the third sound excitation and of the 2D roton, with completion of the first layer, may 
lead to unexpected behaviour of the heat capacity of submonolayer films 



In Fig. |I8|, we show the dynamic structure function associated with the spectra of Fig. [17. 
The strength distribution of the thinner film is more fragmented, as an effect of the finite size 
of the system in the 2;-direction. The relative strength of the ripplon mode is also higher for 
the thinner film, since the surface-to-volume ratio is larger. One notes also that the surface 
mode mantains its identity well above the roton energy. The thicker film has a spectrum 
quite similar to the one of the free surface, apart from the 2D roton at the substrate-helium 
interface, which is still present but affects the spectrum only in the region g ~ 2 A~^. The 
spectrum is dominated by the bulk phonon-roton mode and the lowest surface mode, a non 
negligible strength being also associated with the excited surface modes which lie in between. 

It is worth noticing that the potential between one ^He atom and the H2 substrate is 
similar to that between one ^He and a Graphite substrate coated by two solid ^He layers 
53| , |58[| . Hence, the liquid "^He film on Graphite has a spatial structure similar to the profile 



of Fig. |T0|. Thus we can can safely compare our results with the experimental ones, obtained 
at ILL-Grenoble [p5l, as well as with the variational theory of Ref. |53]. The spectra in 



Fig. are indeed similar to the experimental ones (see for example Fig. 14 of Ref. [Q). 
The energy of the excitations are in rather good agreement, apart from the energy of the 
2D-roton which, however, depends crucially on the substrate-helium potential. A reasonable 
agreement is also obtained for the strength distribution. 
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VI. CONCLUSIONS 



We have presented a new density functional theory for hquid ^He at zero temperature. 
The theory corresponds to an improved version of the density functional introduced in 
Ref. [|I|. It is a phenomenological theory, where known properties of the uniform liquid are 
used to fix the parameters of the functional. The theory is suitable to study inhomogeneous 
states of liquid helium in different geometries like the free surface, droplets, films, bubbles, 
vortices. Equilibrium configurations and excited states can be studied in a unified frame- 
work. Here we have presented results for static and dynamic properties of the free surface 
and films, as well as for the ground state of helium droplets. The comparison with available 
experimental data, as well as with Monte Carlo simulations, is in general satisfactory. The 
present density functional theory turns out to be quantitative even on the scale of inter- 
atomic distances. The improvements with respect to previous density functionals are clearly 
understood in term of the key ingredients of the theory, namely : i) the bulk equation of 
state, ii) the static response function and iii) the phonon-roton dispersion in the uniform 
liquid. In particular we have discussed in details the importance of reproducing the peak of 
the static response function at the roton wavelength in determining static density oscillations 
near the surface and shell structures in helium droplets ("soft packing"). We have also dis- 
cussed the contribution of the collective (one-phonon) modes to sum rules (compressibility 
sum rule, f-sum rule, static form factor), in order to clarify the grounds of TDDF theory, as 
well as to include backfiow effects. The results for the excited states of the free surface and 
films compare quantitatively with experiments. The applications of the same theory to the 
dynamics of droplets and to the phenomenon of quantum evaporation are work in progress. 
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FIGURES 

FIG. 1. Static response function in liquid ^He at zero pressure. Points: experimental data 
|; dotted line: from functional of Refs. jl^Jl^; dashed line: Orsay-Paris functional solid line: 
present functional (Eq. (|0)) 



FIG. 2. Pressure versus density in bulk. Points: experimental data 1 24 1; dashed line: Quantum 



Monte Carlo results |22]; solid line: present density functional. 



FIG. 3. Sound velocity in bulk. Points: experimental data [23j; dashed line: Quantum Monte 



Carlo results |22]; solid line: present density functional. 



FIG. 4. Free surface profile of liquid ^He at zero temperature. Dotted line: functional of Ref. 
1 15, 16 1; dashed line: Orsay-Paris functional [0; solid line: present functional. 



FIG. 5. Density profile of small ^He droplets (density normalized to the bulk value) 
FIG. 6. Density profile of a droplet with 70 ^He atoms. Solid line: present work; triangle: 



DMC simulations of Ref. |^3|; dashed line: Variational (HNC) calculations [34|; circles: DMC 



simulation of Ref. |35] 



FIG. 7. Energy per particle versus N. Solid line: present work; dashed line: results with 



functional of Ref. [||; dotted line: Ref. [16|; crosses: Ref. [37|; circles: Ref. |33| 



FIG. 8. Evaporation energy. Solid line: density functional; dashed line: liquid drop formula. 
Circles: deviation from the liquid drop formula (axis on the right). 

FIG. 9. Evolution of /i vs ^He coverage on several substrates. Solid lines: present work; dashed 
lines: Orsay-Paris functional. The substrates are: H2 (insert) and, from top to bottom, Cs, Rb, 
K, Na, Li, Mg. 
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FIG. 10. Density profile of helium films on solid hydrogen from 0.02 to 0.6 A'^. Solid lines: 
stable films, dotted lines: unstable or metastable films. The growth is not continuous because of 
the two layering transitions related to the formation of the first two layers. 

FIG. 11. Adsorption isotherm on a graphite substrate with T = 0.639 K. Solid line: present 
work; dashed lined: Orsay-Paris functional j^; dots: experimental results extracted from Ref. |44|. 
One layer corresponds to 0.07812 A^^ (i.e., p^^ where po is the bulk density). 

FIG. 12. Phonon-roton dispersion in bulk at zero pressure. Points: experimental data |49|; 
dotted line: Feynman approximation (^) with the experimental static form factor Sq |50|; dashed 
line: Orsay-Paris functional |l],^; dot-dashed line: from Eq. ( [20| ) with TL.^ given in Eq. p8|); solid 
line: same as before but with Ti^ from Eq. 

FIG. 13. Energy weighted moment of the dynamic structure function. Dashed line: total f-sum 



rule; bars: one-phonon contribution as measured in neutron scattering |2C]; solid line: present work 
(Equation (||)). 

FIG. 14. Dispersion relation of bulk and free surface excitations. The threshold A = 8.7 
K, for roton states with different values of parallel wave vector k, is shown as the horizontal 
line. The present result for the surface mode (lowest solid line reaching A at A; = 1.15 A^^) is 
compared with the experimental data on films ||55[ (squares) and with the hydrodynamic dispersion 
of ripplons (short-dashed line). The bulk phonon-roton branch (upper solid line) is compared with 
the experimental one (circles). The threshold for the emission of atoms into the vacuum is also 
shown (dot-dashed curve). 

FIG. 15. Dynamic structure function (in arbitrary units) for = in a slab 50 A thick. The 
lowest line corresponds to A; = 0.3 A~^; the highest one to A; = 1.9 A^^. The dashed lines are the 



phonon-roton and the ripplon dispersion. The (5-function in the definition (62) is replaced by a 
gaussian of width 0.4 K. 
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FIG. 16. Evolution of third sound speed c^s versus coverage on a H2 substrate. Solid line: 
Orsay-Trento results. Dashed line: Orsay-Paris results. 

FIG. 17. Spectrum of films on a H2 substrate for coverages 0.24 A-2 (a) and 0.48 A'^ (b). 
Note the presence of a low energy excitation around 2 A^^ which is associated with a 2D roton 
mode trapped in the first layer. 

FIG. 18. Dynamic structure function (in arbitrary units) associated with the spectra in Fig. 
17, for coverages 0.24 A-2 (a) and 0.48 A-^ (b) on a H2 substrate. The lowest line corresponds to 



/c = 0.1 A ^, the highest one to A: = 2.0 A ^. The (5-function in the definition (63) is replaced by 
a gaussian of width 0.4 K. 
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TABLES 

TABLE I. Values of the potential parameters C3 and D. 





Cs 


Rb 


K 


Na 


Li 


Mg 


Gr 


C3 (KA3) 


673 


754 


812 


1070 


1360 


1775 


2088 


DiK) 


4.41 


4.99 


6.26 


10.4 


17.1 


32.1 


192.6 



TABLE II. Values of the parameters used in Vj(r), see Eq. (46) 



711 


721 


712 


722 




0:2 


-19.7544 


-0.2395 


12.5616 


0.0312 A-2 


1.023 A-2 


0.14912 A-2 
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